def generate_testing_data():
    import soapy
    import os
    import numpy as np
    import matplotlib.pyplot as plt
    import copy
    from astropy.io import fits

    work_dir = os.getcwd()
    config_file = work_dir + '/configuration_dir/generate_data.yaml'
    data_dir = work_dir + '/testing_data'

    sim = soapy.Sim(config_file)
    sim.aoinit()
    sim.makeIMat()


    def get_config():
        """
        This function is used for generating specific parameters for generating test data.
        return: r0, nature guide star magnitude
        """
        r0 = [0.05, 0.06, 0.07, 0.08, 0.09, 0.10, 0.11, 0.12, 0.13, 0.14, 0.15, 0.16, 0.17, 0.18, 0.19, 0.20]
        star_mag = [8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18]
        return [r0, star_mag]


    # Generate training data
    testing_data_element_length = 1000  # The data length for testing of a single configuration parameter
    r0, star_mag = get_config()

    # generate different r0 values firstly
    save_dir = data_dir + '/different_r0'
    GS_mag = 8        # fix the GS_mag to a good value
    for i in r0:
        element_dir = save_dir + '/r0=' + str(i)
        os.makedirs(element_dir)
        for j in range(testing_data_element_length):

            test_data_element_dir = element_dir + '/' + str(j)
            os.makedirs(test_data_element_dir)

            # save_dir = data_dir + '/' + str(i)
            # os.makedirs(save_dir)

            # r0, star_mag = get_config()

            sim.atmos.r0 = i
            sim.atmos.scrnStrengths = (((sim.atmos.r0 ** (-5. / 3.)) * sim.atmos.config.normScrnStrengths) ** (-3. / 5.))
            sim.atmos.randomScrns()

            # get perfect dm command vector
            sim.wfss[0].config.GSMag = 8
            sim.wfss[0].config.photonNoise = False
            perfect_slopes = sim.wfss[0].frame(sim.atmos.scrns)
            perfect_slopes = copy.deepcopy(perfect_slopes)
            perfect_acts = sim.recon.control_matrix.T.dot(perfect_slopes)
            perfect_detector = sim.wfss[0].detector
            perfect_detector = copy.deepcopy(perfect_detector)

            # setting target guide star magnitude
            # sim.wfss[0].config.GSMag = star_mag
            sim.wfss[0].config.photonNoise = True

            distorted_slopes = sim.wfss[0].frame(sim.atmos.scrns)
            distorted_slopes = copy.deepcopy(distorted_slopes)
            distorted_acts = sim.recon.control_matrix.T.dot(distorted_slopes)
            distorted_detector = sim.wfss[0].detector
            distorted_detector = copy.deepcopy(distorted_detector)

            # test correction performance by perfect acts
            dm_screen = sim.dms[0].dmFrame(perfect_acts)
            phase_before_radius = sim.wfss[0].los.frame(sim.atmos.scrns)
            phase_before_radius = copy.deepcopy(phase_before_radius)
            phase_before_radius = phase_before_radius * sim.mask
            phase_after_nm = sim.wfss[0].los.frame(sim.atmos.scrns, dm_screen[np.newaxis, ...])
            phase_after_nm = copy.deepcopy(phase_after_nm)
            phase_after_radius = phase_after_nm * sim.mask * sim.wfss[0].los.phs2Rad
            image_before = sim.sciCams[0].frame(sim.atmos.scrns)
            image_before = copy.deepcopy(image_before)
            inst_strehl_before = sim.sciCams[0].instStrehl
            inst_strehl_before = copy.deepcopy(inst_strehl_before)
            image_after = sim.sciCams[0].frame(sim.atmos.scrns, dm_screen[np.newaxis, ...])
            inst_strehl_after = sim.sciCams[0].instStrehl

            # save result
            title_list = ['phase_before_radius', 'phase_after_radius', 'image_before', 'image_after']
            var_list = [phase_before_radius, phase_after_radius, image_before, image_after]
            plt.figure()
            for _ in range(4):
                plt.subplot(2, 2, _ + 1)
                if _ == 0:
                    plt.imshow(var_list[_], origin='lower')
                    cbar = plt.colorbar()
                    plt.title(title_list[_])
                elif _ == 1:
                    plt.imshow(var_list[_], origin='lower', vmax=cbar.vmax, vmin=cbar.vmin)
                    plt.colorbar()
                    plt.title(title_list[_])
                elif _ == 2:
                    plt.imshow(var_list[_], origin='lower')
                    plt.colorbar()
                    plt.title(title_list[_] + ' with inst_strehl = ' + str(inst_strehl_before))
                else:
                    plt.imshow(var_list[_], origin='lower')
                    plt.colorbar()
                    plt.title(title_list[_] + ' with inst_strehl = ' + str(inst_strehl_after))
            plt.tight_layout()
            plt.savefig(test_data_element_dir + '/ideal_correction_performance')
            plt.close('all')

            # save essential data for training

            hdu = fits.PrimaryHDU(sim.atmos.scrns)
            # hdu.header['r0'] = r0
            # hdu.header['pixel_scale'] = sim.atmos.pixel_scale
            hdulist = fits.HDUList([hdu])
            hdulist.writeto(test_data_element_dir + '/atmos_scrns.fits')
            hdulist.close()

            hdu = fits.PrimaryHDU(phase_before_radius)
            hdu.header['r0'] = i
            hdu.header['pixel_scale'] = sim.atmos.pixel_scale
            hdulist = fits.HDUList([hdu])
            hdulist.writeto(test_data_element_dir + '/phase_before_radius.fits')
            hdulist.close()

            hdu = fits.PrimaryHDU(perfect_detector)
            hdulist = fits.HDUList([hdu])
            hdulist.writeto(test_data_element_dir + '/perfect_detector.fits')
            hdulist.close()

            hdu = fits.PrimaryHDU(perfect_slopes)
            hdulist = fits.HDUList([hdu])
            hdulist.writeto(test_data_element_dir + '/perfect_slopes.fits')
            hdulist.close()

            hdu = fits.PrimaryHDU(perfect_acts)
            hdulist = fits.HDUList([hdu])
            hdulist.writeto(test_data_element_dir + '/perfect_acts.fits')
            hdulist.close()

            hdu = fits.PrimaryHDU(distorted_detector)
            hdulist = fits.HDUList([hdu])
            hdulist.writeto(test_data_element_dir + '/distorted_detector.fits')
            hdulist.close()

            hdu = fits.PrimaryHDU(distorted_slopes)
            hdulist = fits.HDUList([hdu])
            hdulist.writeto(test_data_element_dir + '/distorted_slopes.fits')
            hdulist.close()

            hdu = fits.PrimaryHDU(distorted_acts)
            hdulist = fits.HDUList([hdu])
            hdulist.writeto(test_data_element_dir + '/distorted_acts.fits')
            hdulist.close()

            print(f'Firstly different r0: r0={i} {j} data has completed.')

            # sim.wfss[0].addPhotonNoise()

    # generate different guide star magnitude data for testing secondly
    save_dir = data_dir + '/different_mag'
    # fix the r0 of the atmosphere to a good value
    sim.atmos.r0 = 0.2
    sim.atmos.scrnStrengths = (((sim.atmos.r0 ** (-5. / 3.)) * sim.atmos.config.normScrnStrengths) ** (-3. / 5.))
    for i in star_mag:
        element_dir = save_dir + '/star_mag=' + str(i)
        os.makedirs(element_dir)
        for j in range(testing_data_element_length):

            test_data_element_dir = element_dir + '/' + str(j)
            os.makedirs(test_data_element_dir)

            # save_dir = data_dir + '/' + str(i)
            # os.makedirs(save_dir)

            # r0, star_mag = get_config()

            # sim.atmos.r0 = i
            # sim.atmos.scrnStrengths = (((sim.atmos.r0 ** (-5. / 3.)) * sim.atmos.config.normScrnStrengths) ** (-3. / 5.))
            sim.atmos.randomScrns()

            # get perfect dm command vector
            sim.wfss[0].config.GSMag = 8
            sim.wfss[0].config.photonNoise = False
            perfect_slopes = sim.wfss[0].frame(sim.atmos.scrns)
            perfect_slopes = copy.deepcopy(perfect_slopes)
            perfect_acts = sim.recon.control_matrix.T.dot(perfect_slopes)
            perfect_detector = sim.wfss[0].detector
            perfect_detector = copy.deepcopy(perfect_detector)

            # setting target guide star magnitude
            sim.wfss[0].config.GSMag = i
            sim.wfss[0].config.photonNoise = True

            distorted_slopes = sim.wfss[0].frame(sim.atmos.scrns)
            distorted_slopes = copy.deepcopy(distorted_slopes)
            distorted_acts = sim.recon.control_matrix.T.dot(distorted_slopes)
            distorted_detector = sim.wfss[0].detector
            distorted_detector = copy.deepcopy(distorted_detector)

            # test correction performance by perfect acts
            dm_screen = sim.dms[0].dmFrame(perfect_acts)
            phase_before_radius = sim.wfss[0].los.frame(sim.atmos.scrns)
            phase_before_radius = copy.deepcopy(phase_before_radius)
            phase_before_radius = phase_before_radius * sim.mask
            phase_after_nm = sim.wfss[0].los.frame(sim.atmos.scrns, dm_screen[np.newaxis, ...])
            phase_after_nm = copy.deepcopy(phase_after_nm)
            phase_after_radius = phase_after_nm * sim.mask * sim.wfss[0].los.phs2Rad
            image_before = sim.sciCams[0].frame(sim.atmos.scrns)
            image_before = copy.deepcopy(image_before)
            inst_strehl_before = sim.sciCams[0].instStrehl
            inst_strehl_before = copy.deepcopy(inst_strehl_before)
            image_after = sim.sciCams[0].frame(sim.atmos.scrns, dm_screen[np.newaxis, ...])
            inst_strehl_after = sim.sciCams[0].instStrehl

            # save result
            title_list = ['phase_before_radius', 'phase_after_radius', 'image_before', 'image_after']
            var_list = [phase_before_radius, phase_after_radius, image_before, image_after]
            plt.figure()
            for _ in range(4):
                plt.subplot(2, 2, _ + 1)
                if _ == 0:
                    plt.imshow(var_list[_], origin='lower')
                    cbar = plt.colorbar()
                    plt.title(title_list[_])
                elif _ == 1:
                    plt.imshow(var_list[_], origin='lower', vmax=cbar.vmax, vmin=cbar.vmin)
                    plt.colorbar()
                    plt.title(title_list[_])
                elif _ == 2:
                    plt.imshow(var_list[_], origin='lower')
                    plt.colorbar()
                    plt.title(title_list[_] + ' with inst_strehl = ' + str(inst_strehl_before))
                else:
                    plt.imshow(var_list[_], origin='lower')
                    plt.colorbar()
                    plt.title(title_list[_] + ' with inst_strehl = ' + str(inst_strehl_after))
            plt.tight_layout()
            plt.savefig(test_data_element_dir + '/ideal_correction_performance')
            plt.close('all')

            # save essential data for training

            hdu = fits.PrimaryHDU(sim.atmos.scrns)
            # hdu.header['r0'] = r0
            # hdu.header['pixel_scale'] = sim.atmos.pixel_scale
            hdulist = fits.HDUList([hdu])
            hdulist.writeto(test_data_element_dir + '/atmos_scrns.fits')
            hdulist.close()

            hdu = fits.PrimaryHDU(phase_before_radius)
            hdu.header['r0'] = sim.atmos.r0
            hdu.header['pixel_scale'] = sim.atmos.pixel_scale
            hdulist = fits.HDUList([hdu])
            hdulist.writeto(test_data_element_dir + '/phase_before_radius.fits')
            hdulist.close()

            hdu = fits.PrimaryHDU(perfect_detector)
            hdulist = fits.HDUList([hdu])
            hdulist.writeto(test_data_element_dir + '/perfect_detector.fits')
            hdulist.close()

            hdu = fits.PrimaryHDU(perfect_slopes)
            hdulist = fits.HDUList([hdu])
            hdulist.writeto(test_data_element_dir + '/perfect_slopes.fits')
            hdulist.close()

            hdu = fits.PrimaryHDU(perfect_acts)
            hdulist = fits.HDUList([hdu])
            hdulist.writeto(test_data_element_dir + '/perfect_acts.fits')
            hdulist.close()

            hdu = fits.PrimaryHDU(distorted_detector)
            hdulist = fits.HDUList([hdu])
            hdulist.writeto(test_data_element_dir + '/distorted_detector.fits')
            hdulist.close()

            hdu = fits.PrimaryHDU(distorted_slopes)
            hdulist = fits.HDUList([hdu])
            hdulist.writeto(test_data_element_dir + '/distorted_slopes.fits')
            hdulist.close()

            hdu = fits.PrimaryHDU(distorted_acts)
            hdulist = fits.HDUList([hdu])
            hdulist.writeto(test_data_element_dir + '/distorted_acts.fits')
            hdulist.close()

            print(f'Secondly different mag: mag={i} {j} data has completed.')

            # sim.wfss[0].addPhotonNoise()

    return True
